function [Yield,Iterations] = bndyield(Price, CouponRate, Settle, Maturity,varargin)
%BNDYIELD Yield to maturity for a fixed income security.
%   Given NumBonds bonds with date parameters and clean prices, this
%   function returns the yields to maturity.
%
%   Yield = bndyield(Price, CouponRate, Settle, Maturity)
%
%   Yield = bndyield(Price, CouponRate, Settle, Maturity, Period,
%      Basis, EndMonthRule, IssueDate, FirstCouponDate, LastCouponDate,
%      StartDate, Face)
%
%   Yield = bndyield(Price, CouponRate, Settle, Maturity,
%      'Param1','Value1',...)
%
%   Optional Inputs: Period, Basis, EndMonthRule, IssueDate, FirstCouponDate,
%                    LastCouponDate, StartDate, Face, LastCouponInterest,
%                    CompoundingFrequency, DiscountBasis
%
%   Note:
%   - For SIA conventions, the Price and Yield are related by the formula:
%     Price + Accrued Interest = sum(Cash_Flow*(1+Yield/2)^(-Time))
%     where the sum is over the bond's cash flows and corresponding times in
%     units of semi-annual coupon periods.
%     For ISMA conventions, the Price and Yield are related by the formula:
%     Price + Accrued Interest = sum(Cash_Flow*(1+Yield)^(-Time))
%   - All non-scalar or empty matrix input arguments must be either NUMBONDSx1
%     or 1xNUMBONDS conforming vectors.
%   - Fill unspecified entries in input vectors with NaN.
%   - Dates can be serial date numbers or date strings.
%   - Optional inputs can be specified as parameter value pairs.  If
%     LastCouponInterest, CompoundingFrequency or DiscountBasis are input,
%     optional inputs must be specified as parameter value pairs.
%     Otherwise, optional inputs may be specified by order according to the
%     help.
%
%   Inputs:
%        Price - Clean price.
%
%   CouponRate - Coupon rate in decimal form.
%
%       Settle - Settlement date.
%
%     Maturity - Maturity date.
%
%   Optional Inputs:
%            Period - Number of coupons payments per year.
%                     Possible values include:
%                     0, 1, 2 (default), 3, 4, 6, 12
%
%             Basis - Day-count basis.
%                     Possible values include:
%                     0 - actual/actual (default)
%                     1 - 30/360 SIA
%                     2 - actual/360
%                     3 - actual/365
%                     4 - 30/360 PSA
%                     5 - 30/360 ISDA
%                     6 - 30/360 European
%                     7 - actual/365 Japanese
%                     8 - actual/actual ISMA
%                     9 - actual/360 ISMA
%                    10 - actual/365 ISMA
%                    11 - 30/360 ISMA
%                    12 - actual/365 ISDA
%                    13 - bus/252
%
%      EndMonthRule - End-of-month rule; default is 1 (in effect)
%                     0 - Rule is NOT in effect for the bond(s)
%                     1 - (default) Rule is in effect for the bond(s) (meaning
%                         that a security that pays coupon interest on the last
%                         day of the month will always make payment on the last
%                         day of the month)
%
%         IssueDate - Bond issue date.
%
%   FirstCouponDate - Irregular or normal first coupon date.
%
%    LastCouponDate - Irregular or normal last coupon date.
%
%         StartDate - Forward starting date of payments.
%
%              Face - Face value of the bond; default is 100.
%
%   LastCouponInterest - Compounding convention for computing the yield of
%                        a bond in the last coupon period, i.e.: with only
%                        the last coupon and the face value to be repaid.
%                        Choices are 'simple' or 'compound'.
%
%   CompoundingFrequency - Compounding frequency for yield calculation.  By
%                          default, SIA bases (0-7) and BUS/252 use a semi-annual
%                          compounding convention and ISMA bases (8-12) use
%                          an annual compounding convention.
%
%   DiscountBasis - Basis used to compute the discount factors for
%                   computing the yield.  The default behavior is for SIA
%                   bases to use the actual/actual day count to compute
%                   discount factors, and for ISMA day counts and BUS/252
%                   to use the specified basis.
%
%   Outputs:
%   Yield - [NUMBONDSx1 vector] of the yield to maturity with semi-annual
%           compounding.
%
%   See also BNDPRICE, CFAMOUNTS.

%   Copyright 1995-2009 The MathWorks, Inc.
%   $Revision: 1.13.2.15 $   $Date: 2009/11/05 16:59:02 $

% Checking input arguments
if nargin < 4
    error('finance:bndyield:notEnoughInputs', 'Too few inputs.');
end

% Check to see whether we have the case of ordered inputs or PV pairs
if ~isempty(varargin)
    if ischar(varargin{1})
        
        p = inputParser;
        
        p.addParamValue('discountbasis',NaN,@isvalidbasis);
        p.addParamValue('compoundingfrequency',NaN,@(x) ismember(x,[-1 1 2 3 4 6 12]));
        p.addParamValue('lastcouponinterest',{'compound'},...
            @(x) ismember(x,{'compound','simple'}));
        
        p.addParamValue('period',[]);
        p.addParamValue('basis',[]);
        p.addParamValue('endmonthrule',[]);
        p.addParamValue('issuedate', []);
        p.addParamValue('firstcoupondate', []);
        p.addParamValue('lastcoupondate', []);
        p.addParamValue('startdate', []);
        p.addParamValue('face',[]);
        
        try
            p.parse(varargin{:});
        catch ME
            newME = MException('finance:bndyield:optionalInputError',...
                'Error in input arguments');
            newME = addCause(newME,ME);
            throw(newME)
        end
        
        Frequency = p.Results.compoundingfrequency;
        LastCouponInterest = cellstr(p.Results.lastcouponinterest);
        Basis =  p.Results.basis;
        DiscountBasis = p.Results.discountbasis;
        
        try
            [CouponRate, Settle, Maturity,~,Basis] = instargbond(CouponRate,...
                Settle, Maturity,[], Basis);
        catch ME
            throwAsCaller(ME)
        end
    else
        Frequency = NaN;
        LastCouponInterest = {'compound'};
        DiscountBasis = NaN;
        
        try
            [CouponRate, Settle, Maturity,~,Basis] = instargbond(CouponRate,...
                Settle, Maturity,varargin{:});
        catch ME
            throwAsCaller(ME)
        end
    end
else
    
    Frequency = NaN;
    LastCouponInterest = {'compound'};
    DiscountBasis = NaN;
    
    try
        [CouponRate, Settle, Maturity,~,Basis] = instargbond(CouponRate, Settle, Maturity);
    catch ME
        throwAsCaller(ME)
    end
end

% Strip out lastcouponinterest
if any(strcmpi(varargin,'lastcouponinterest'))
    lciidx = find(strcmpi(varargin,'lastcouponinterest'));
    varargin([lciidx lciidx+1]) = [];
end

% Scalar expansion on non-bond inputs
try
    [CouponRate, Settle, Maturity, Price, Frequency, LastCouponInterest, DiscountBasis] = ...
        finargsz(1, CouponRate, Settle, Maturity, Price, Frequency, ...
        LastCouponInterest, DiscountBasis);
catch ME
    throwAsCaller(ME)
end
% Generate cash flows
[CFlowAmounts, ~, TFactors] = cfamounts(CouponRate, Settle, Maturity, varargin{:});

nBonds = length(Price);

% If no compounding frequency specified, use 2 for SIA bases, 1 for ISMA.
if isnan(Frequency)
    Frequency = 2*ones(nBonds,1);
    if isnan(DiscountBasis)
        i = isisma(Basis);
    else
        i = isisma(DiscountBasis);
    end
    Frequency(i) = 1;
end

CFlowAmounts(:, 1) = CFlowAmounts(:, 1) - Price;

% Initial Guess for the discount based on Yield = Coupon Rate
CouponRate(isnan(CouponRate)) = 0;
XGuess = 1./(1 + CouponRate./Frequency);

fZeroOptions = optimset('fzero');
fZeroOptions = optimset(fZeroOptions, 'TolX', 1e-12,'Display', 'off');

fSolveOptions = optimset('fsolve');
fSolveOptions = optimset(fSolveOptions, 'TolX', 1e-12,'TolFun', 1e-12,...
    'Display', 'off','LargeScale', 'off');

Yield = zeros(nBonds,1);
Iterations = NaN*ones(nBonds,1);

% Calculate yield
for i = 1:nBonds
    if isnan(Price(i))
        Yield(i) = NaN;
    elseif sum(~isnan(CFlowAmounts(i,:))) == 2 % Settlement in last coupon period
        % Solve directly for PerDisc
        switch lower(LastCouponInterest{i})
            case 'compound'
                PerDisc = (-CFlowAmounts(i,1)/CFlowAmounts(i,2))^(1/TFactors(i,2));
                Yield(i) = (1./PerDisc - 1).*Frequency(i);
                Flag = 1;
            case 'simple'
                Yield(i) = (-CFlowAmounts(i,2)./CFlowAmounts(i,1)-1)./TFactors(i,2)...
                    .*Frequency(i);
                Flag = 1;
        end
    else % General case, settlement before the last coupon period
        
        objfun = @(x) nansum(CFlowAmounts(i,:) .* (x.^TFactors(i,:)));
        [PerDisc,~, Flag, Output] = fzero(objfun,XGuess(i), fZeroOptions);
        
        % if fzero was unable to find a zero, try using fsolve
        if Flag < 0
            [PerDisc, ~, Flag, Output] = fsolve(objfun, XGuess(i), fSolveOptions);
        end
        
        Yield(i) = (1./PerDisc - 1).*Frequency(i);
        Iterations(i) = Output.iterations-1;
    end
end % end of loop over bonds

if Flag < 0
    warning('finance:bndyield:solutionConvergenceFailure', ...
        'Could not solve for the yield.\n');
end